Delaunay Triangulations in Linear Time ? ( Part I ) ∗ Kevin Buchin

We present a new and simple randomized algorithm for constructing the Delaunay triangulation using nearest neighbor graphs for point location. It runs in linear expected time for points in the plane with polynomially bounded spread, i.e., if the ratio between the largest and smallest pointwise distance is polynomially bounded. This also holds for point sets with bounded spread in higher dimensions as long as the expected complexity of the Delaunay triangulation of a sample of the points is linear in the sample size. Chan and P atra3cu [6, 7] presented o(N log N) randomized algorithms for constructing Voronoi Diagrams of points in the plane (from which the Delaunay triangulation can be computed in linear time and vice-versa) under suitable models of computation. We improve on these results by presenting an O(N) randomized algorithm for the Delaunay triangulation in the plane. The algorithm is not restricted to two dimensions and it runs in linear expected time as long as the expected complexity of the Delaunay triangulation of a random sample of the input points is linear in the sample size. An example of linear complexity Delaunay triangulation are suitably sampled (d− 1)-dimensional polyhedra in IR (but no guarantees are given for samples). Our algorithm locates points by combining the history (i.e., the Delaunay tree [2, 3]) of a randomized incremental construction with a sequence of nearest neighbor graph computations. For the nearest neighbor graphs we use a recent result by Chan [5] that links well-separated pair decompositions to sorting. By the use of radix sort this results in a linear time algorithm for well-separated pair decompositions and as a consequence for nearest neighbor graphs. We will use the same assumptions as Chan on the model of computation and the point set. The model of computation is a real-RAM with the oor function available and a word size of at least log N . The input point set should have polynomially bounded spread, i.e., the ratio of the largest and smallest point to point distance should be bounded by a polynomial in the size of the point set. But also other combinations of models of computation and sorting algorithms can be used. ∗This is the rst part of a longer paper that additionally contains results from [4] This research was supported by the Deutsche Forschungsgemeinschaft within the European graduate program 'Combinatorics, Geometry, and Computation' (No. GRK 588/2) and by the Netherlands' Organisation for Scienti c Research (NWO) under BRICKS/FOCUS project no. 642.065.503. †Dept. of Information and Comp. Sci., Utrecht Univ., Netherlands; [email protected]. 1 ar X iv :0 81 2. 03 87 v1 [ cs .C G ] 1 D ec 2 00 8 Algorithm 1: Incremental Construction using Nearest Neighbor Graph Input: Finite point set P in IR Output: Delaunay triangulation of P 1 Split P into rounds R1, . . . Rm of doubling size with R1 of constant size and set Sj ← ⋃ 1≤i≤j Ri (j = 1 . . . m). 2 Insert points in R1 into the Delaunay triangulation using history for point location. 3 For k = 2, . . . ,m insert points in Rk into the Delaunay triangulation as follows: 3.1 Set Tk−1 ← Rk, Ti ← ∅ (0 ≤ i < k − 1), and j ← k − 1. 3.2 While Tj 6= ∅ and j > 0: compute NNG(Tj ∪ Sj) and from each connected component with no vertex in Sj add the rst point of the component to Tj−1; then set j ← j − 1. 3.3 While j < k − 1: locate Tj (if not empty) in DT(Sj+1) using history starting at DT(Sj); then locate Tj+1 in DT(Sj+1) by walking through the connected components starting at an already located point; then set j ← j + 1 General Setup. We construct the Delaunay triangulation of a nite point set P ⊂ IR by a randomized incremental construction using a history. The point location is accelerated by locating points at intermediate levels in the history instead of the top, see Algorithm 1. Given an insertion order we group the points into rounds R1, . . . , Rm in accordance with the order, i.e., the points in Ri are in the insertion order before the points in Ri+1 for 1 ≤ i < m. The rounds double in size, i.e., |R1| is constant, and |Ri+1| = 2|Ri| (with possibly the exception of the last round for which |Rm| ≤ 2Rm−1). Let Sj := ⋃ 1≤i≤j Ri denote the points inserted in or before round j. Additionally to the history we store the Delaunay triangulations of the Sj. Note that the rounds are only used for facilitating the point location; the insertion order remains the same. Point Location in a Round. The points of the rst round are located in the standard way using the history. At the beginning of round k (2 ≤ k ≤ m) the points of the round Rk are located in the Delaunay triangulation of Sk−1 in the following way: Let Tk−1 := Rk. We compute the nearest neighbor graph of Tk−1 ∪ Sk−1. For connected components of the nearest neighbor graph without a vertex in Sk−1 we include the rst point (according to the insertion order) of the component in a set Tk−2. We repeat the same procedure higher up in the history, i.e., we compute the nearest neighbor graph of Tk−2 ∪ Sk−2, for each connected component without a vertex in Sk−2 we include the rst point in a set Tk−3, and so on. We stop this process with the construction of T0 (or earlier with Tj if Tj−1 is empty. For simplicity we describe the algorithm for the case that T0 is not empty). This yields a hierarchy of sets Tk−1 ⊃ Tk−2 ⊃ · · · ⊃ T0. Now we locate the points in T0 in DT(S1) by using the history, i.e., we use the history to nd a con icting simplex and then locally search for the simplex containing T0. For locating T1 we have the following situation: each connected component of the nearest neighbor graph of T1 ∪ S1 either has a vertex in S1 or has a vertex in T0, thus each component has a vertex already located in DT(S1). We traverse each component starting at an already located vertex,